%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%
% Fonction de résolution du problème flux-free sur la poutre
%
% Entr�es :
%
% sous-maillage : vecteur de coordonnées des noeuds du sous-maillage
% uh : [vecteur 4x1] extraction de la solution éléments finis
% epsilon : =1 sur le premier �l�ment du patch,
%           =-1 sur le deuxi�me �l�ment du patch,
%
% Sorties :
%
% ie : [matrice 4x4] matrice ie

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

function [ normeErreur,erreur,F ] = Resolution_fluxfree( donnees , nSousElements , maillage , uh )


sousPoutre.L=donnees.poutre.L/nSousElements;

erreur = zeros(2*nSousElements*(size(maillage,1)-1)+2,1) ;
F = zeros(2*nSousElements*(size(maillage,1)-1)+2,1) ;


% Création de la matrice de sousmaillage

sousMaillage = zeros( (size(maillage,1)-1)*nSousElements+1 , 1);

for i=1:size(maillage,1)
    sousMaillage( (i-1)*nSousElements+1 ) = maillage(i) ;
end

for i=1:size(maillage,1)-1
    
    coeff = ( sousMaillage(i*nSousElements+1 , 1 ) -  sousMaillage( (i-1)*nSousElements+1 , 1))/nSousElements ;
    for k=1:nSousElements-1
        sousMaillage((i-1)*nSousElements+1+k,1)=sousMaillage((i-1)*nSousElements+1,1)+coeff*k ;
    end
    
end

%%%%%%%%




% Calcul de la matrice de rigidité de l'ensemble (+ CL)

K = AssemblageRigidite(sousMaillage,donnees.poutre) ;


C = [ eye(2), zeros(2,size(K,1)-2) ];
Kcl = [K , C' ; C , zeros(2) ] ;

%%%%%%%%%



for patch=1:donnees.nombreElements+1
        

    
    % Calcul de Fe pour le patch
    
    Fcl = Fe_patch_ff(patch,donnees.nombreElements+1,donnees,nSousElements,maillage,uh,donnees.poutre,sousPoutre) ;
    
%     plot(Fcl)
%     pause
    
    %%%%%%%%%%%%
    
    
    % Résolution du problème pour un patch
    
    Eecl = Kcl\Fcl ;
    
    Ee = Eecl(1:size(Eecl,1)-2) ;
%     plot(Ee)
%     pause
    
%     [Ee,Fe] = Resolution_patch_fluxfree(patch,donnees.nombreElements,donnees,nSousElements,maillage,uh,donnees.poutre,sousPoutre) ;
    
    %%%%%%%%%
    
    
%     F=F+Fcl(1:size(Fcl,1)-2);
    
    erreur = erreur + Ee ;
    
end
   
    normeErreur = erreur'*K*erreur ;
    
%  
%     
%     if patch==1
%         
%         erreur(1:(nSousElements+1)*2) =  Ee ;
% %        F(1:(nSousElements+1)*2) =  Fe ;
%         
%     elseif patch==donnees.nombreElements
%         
%         erreur(size(U,1)-(nSousElements+1)*2+1 : size(U,1)) =  Ee ;
% %        F(size(F,1)-(nSousElements+1)*2+1 : size(F,1)) =  Fe ;
%         
%     else
%         
%         erreur( (2*nSousElements-1)*(patch-2)-(patch-3) : (2*nSousElements-1)*(patch-2)-(patch-3) + size(Ue,1)-1 ) = erreur( (2*nSousElements-1)*(patch-2)-(patch-3) : (2*nSousElements-1)*(patch-2)-(patch-3) + size(Ue,1)-1 ) +  Ee ;
% %        F( (2*nSousElements-1)*(patch-2)-(patch-3) : (2*nSousElements-1)*(patch-2)-(patch-3) + size(Fe,1)-1 ) = F( (2*nSousElements-1)*(patch-2)-(patch-3) : (2*nSousElements-1)*(patch-2)-(patch-3) + size(Fe,1)-1 ) + Fe ;
%         
%         
%     end
%     
    
    
    
    
    
    
    
    %     patch
    %
    %     Ue
    %     Fe
    %
    %     pause
    %
%     normeU = F'*U ;
    
    
    
end